Tectocoris analysis
Kei-Lin Ooi
In this study, we examined the role of reflectivity in thermoregulation and its influence on warning signal variation in the cotton harlequin bug. Here, we have provided the code used to inspect and extract the results of the reflectivity measurements and heating experiments, as well as the code used to analyse the relationships between reflectivity, body size, and heating.
Setup
Load our packages!
library(pavo)
library(tidyverse)
library(corrplot)
library(readxl)
library(writexl)
library(lattice)
library(modelsummary)
library(tinytable)
library(ggimage)Retrieve data
Retrieve spectrophotometry, irradiance, and filter transmittance data.
### --- Spectrophotometry data from bugs --- ###
specs <- read_excel("data/specs.xlsx", sheet = "bug_specs") %>% as.rspec()
# specs <- read.csv("Dryad/specs.csv") %>% as.rspec()
#### --- Irradiance data --- ####
# Standard Solar Irradiance from 280 to 4000 nm
# ASTM G173-03 Reference Spectra Derived from SMARTS v. 2.9.2 (AM1.5)
# Global tilt W*m-2*nm-1
# from https://www.pveducation.org/pvcdrom/appendices/standard-solar-spectra
sun <- read_excel("data/specs.xlsx", sheet = "irradiance") %>%
as.rspec(., whichwl = "wl") %>%
filter(wl>=400 & wl<=1700) # Limit irradiance to the appropriate wavelength range for experiment
#### --- Filter tranmittance data --- ####
filters<-read_excel("data/specs.xlsx", sheet = "filter_transmittance") %>%
as.rspec(.) %>%
procspec(.,fixneg = "zero") %>%
rename(Full.f = halfsun) %>% # rename
mutate(vis.f = Full.f*visfilterTransmittance/100) %>% # vis filter
mutate(nir.f = Full.f*nirfilterTransmittance/100) %>% # nir filter
select(-visfilterTransmittance) %>%
select(-nirfilterTransmittance)Retrieve bug and heating data.
# get reference datasheet
datasheet <- read_excel("data/tectocoris_data.xlsx")
# get vector of individual bugs
bugs <- dir(path = "data/heating_data", pattern = "CSV") %>% gsub(pattern="\\w+\\_(\\w+).CSV", replacement="\\1") %>% unique()
# get treatments
treatments <- c("full", "nir", "uvvis")Reflectivity analysis
We measured the total hemispherical reflectance for a subset of individuals (n=11) in order to calculate the reflectivity of each bug. Looking at the effect of reflectivity on heating instead of colour allows us to investigate the effect of full spectrum solar light.
Compiling spectrometer data
When using two spectrometers, sometimes there is a jump/mismatch where they join, often because the specimen being measured is not flat and there light bleeding in through the probe. Here is code to correct.
#### --- Correct jump in reflectance curve --- ####
# new dataframe to save the fixed specs into
specs2 <- data.frame(
wl = specs$wl)
# loop to fix all specs
for(i in 2:ncol(specs)){
# calculate the difference between each reflectance measurement and the previous measurement
lag1 <- specs[,i] - lag(specs[,i])
# create a vector that identifies rows where the difference is greater than 1 or less than -1 (i.e. the jump)
lagcondition <- lag1 > 1 | lag1 < -1
# identifies the position in the dataset, so we can cut at a particular wavelengths
positions <- which(lagcondition == TRUE)
positions <- positions[positions >= 100] # set it to not fix jumps at the early wavelengths (<500)
# *** THIS IS A CHECK! LOOK AT THIS OUTPUT *** it shows which wavelengths will be cut
cutwls <- specs$wl[positions]
print(colnames(specs)[i])
print(cutwls)
if (length(cutwls) == 0) {
specs2[,i] <- specs[,i]
} else {
# calculate difference - takes an average before and after the drop, in case there is noise
loweravg <- seq(min(positions)-5, min(positions)-1)
upperavg <- seq(max(positions)+1, max(positions)+5)
total <- mean(specs[loweravg, i]) - mean(specs[upperavg, i])
# final adjustment - visible and UV wavelengths are not altered, but anything above the jump is altered
# the wavelengths where the jump occurred are removed
specs2[,i] <- case_when(
specs$wl < specs$wl[min(positions)] ~ specs[,i],
specs$wl > specs$wl[max(positions)] ~ specs[,i] + total,
specs$wl == specs$wl[positions] ~ NA_integer_)
}
}
# replace the column names
colnames(specs2) <- colnames(specs)
# turn back to rspec for easy plotting to check the final result
specs2 <- as.rspec(specs2)
#### --- Smooth curves --- ####
# # this will produce graphs with 5 levels of smoothing. Pick the best, and input the value into 'span' below.
# plotsmooth(specs2, minsmooth = 0.05, maxsmooth = 0.5, curves= 5, ask = FALSE)
specs3 <- procspec(specs2, opt='smooth', fixneg="addmin", span=0.05)
#### --- Pivot data into format appropriate for ggplot --- ####
specs_longer <- specs3 %>%
# Make dataset long
pivot_longer(2:length(specs), names_to = "specID", values_to = "Refl") %>% # change the number to how many files u have
# Create new columns for each variable based on spec name
mutate(code = str_split(specID, pattern = "_", simplify = TRUE)[,1],
region_colour = str_split(specID, pattern = "_", simplify = TRUE)[,3],
context = str_split(specID, pattern = "_", simplify = TRUE)[,2],
replicate = str_split(specID, pattern = "_", simplify = TRUE)[,6])
# add colour column, considering sex
specs_longer$code_replicate <- paste0(specs_longer$code, "_", specs_longer$replicate)
specs_longer$sex <- gsub("[0-9]", "", specs_longer$code)
specs_longer$sex_colour <- paste0(specs_longer$sex, "_", specs_longer$region_colour)
specs_longer$sex <- factor(specs_longer$sex,
levels = c("male", "female"))Visualise spectrometer readings.
# plot spec curves
plot(specs) # original
plot(specs2) # fixed jump
plot(specs3) # smoothed
# plot reflectance reads by individual
plotcolours <- c("black"="black", "blue"="blue", "green"="darkgreen", "orange" = "orange", "red" = "red")
ggplot(specs_longer,
# %>% filter(context == "dead"),
aes(x = wl, y = Refl)) +
geom_line(aes(color = region_colour), alpha = 0.3) +
facet_wrap(~ code, nrow = 3) +
# coord_cartesian(ylim = c(0, 110), xlim = c(900, 1000)) +
scale_colour_manual(values = plotcolours) +
theme_bw()
Averaging reflectance
Average all reads of a colour for each individual bug.
# define the function for calculating the standard error of the mean
se <- function(x) {
sd(x) / sqrt(length(x))
}
# average all reads of a colour for each individual bug
specs_avg <- specs_longer %>%
group_by(wl, code, sex_colour) %>%
summarise(mean_refl = mean(Refl))Calculating reflectivity
1. Data frame
First, create a data frame with reflectance and irradiance, one per colour and for each treatment (full, NIR, and visible).
We need the reflectivity for the iridescent portion, as well as for the non-iridescent portion for males and females. We average blue + green to represent iridescence, and red + orange to represent non-iridescence. Males and females non-iridescent patches are averaged separately because males are more red and females are more orange. Note that we ignore the black colour category since it is a rarer colour that some males exhibit in their iridescent patches.
# get average across colour categories (blue and green averaged, red and orange averaged separately for males/females)
full_iri <- specs_avg %>%
filter(sex_colour=="male_blue" | sex_colour=="male_green" | sex_colour=="female_blue") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl), # calculate the mean reflectance
se_grouped_refl = se(mean_refl)) %>% # calculate standard error of the mean
data.frame(filters$Full.f) %>% # add the irradiance values
mutate(colour = "iridescent") # add a column indicating colour
full_noniri_male <- specs_avg %>%
filter(sex_colour=="male_red" | sex_colour=="male_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$Full.f) %>%
mutate(colour = "male non-iridescent")
full_noniri_female <- specs_avg %>%
filter(sex_colour=="female_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$Full.f) %>%
mutate(colour = "female non-iridescent")
# for nir
nir_iri <- specs_avg %>%
filter(sex_colour=="male_blue" | sex_colour=="male_green" | sex_colour=="female_blue") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl), # calculate the mean reflectance
se_grouped_refl = se(mean_refl)) %>% # calculate standard error of the mean
data.frame(filters$nir.f) %>% # add the irradiance values
mutate(colour = "iridescent") # add a column indicating colour
nir_noniri_male <- specs_avg %>%
filter(sex_colour=="male_red" | sex_colour=="male_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$nir.f) %>%
mutate(colour = "male non-iridescent")
nir_noniri_female <- specs_avg %>%
filter(sex_colour=="female_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$nir.f) %>%
mutate(colour = "female non-iridescent")
# for visible
vis_iri <- specs_avg %>%
filter(sex_colour=="male_blue" | sex_colour=="male_green" | sex_colour=="female_blue") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl), # calculate the mean reflectance
se_grouped_refl = se(mean_refl)) %>% # calculate standard error of the mean
data.frame(filters$vis.f) %>% # add the irradiance values
mutate(colour = "iridescent") # add a column indicating colour
vis_noniri_male <- specs_avg %>%
filter(sex_colour=="male_red" | sex_colour=="male_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$vis.f) %>%
mutate(colour = "male non-iridescent")
vis_noniri_female <- specs_avg %>%
filter(sex_colour=="female_orange") %>%
group_by(wl) %>%
summarise(
mean_grouped_refl = mean(mean_refl),
se_grouped_refl = se(mean_refl)) %>%
data.frame(filters$vis.f) %>%
mutate(colour = "female non-iridescent")Plot mean and SE for iridescent and non-iridescent patches.
full_group <- rbind(full_iri, full_noniri_male, full_noniri_female)
# Plot mean and SE for iridescent and non-iridescent patches
full_group$colour <- factor(full_group$colour,
levels = c("female non-iridescent","male non-iridescent", "iridescent"))
# plot the means and the standard error of the mean
ggplot(full_group, aes(x = wl, y = mean_grouped_refl, colour = colour)) +
geom_line(alpha = 0.7, linewidth = 1) + # Line with colour mapped globally
geom_ribbon(aes(ymin = mean_grouped_refl - se_grouped_refl, # add the standard error of the mean
ymax = mean_grouped_refl + se_grouped_refl, fill = colour),
alpha = .2, colour = NA) + # Fill for ribbon, no border line
geom_vline(xintercept = 700, linetype = 2, alpha = .5) + # add the border for visual spectra
annotate("text", x = 700, y = Inf, label = "near infrared →",
vjust=1.35, hjust=-0.1, alpha=0.5) + # Label for NIR line
labs(x = "Wavelength", y = "Reflectance",
colour = "Colour", fill = "Colour") + # Same label for both fill and colour
scale_colour_manual(values = c("iridescent" = "blue",
"female non-iridescent" = "orange",
"male non-iridescent" = "red")) +
scale_fill_manual(values = c("iridescent" = "blue",
"female non-iridescent" = "orange",
"male non-iridescent" = "red")) +
theme_bw()
2. Multiply
Find the product between reflectivity and irradiance. Note: Here, reflectance is divided by 100, i.e. as a fraction that can vary between 0 and 1. We do this to be able to illustrate the results in a plot. However, it is not required to divide by 100. If we divide by 100, the answer will be the reflectivity as a proportion, and if we do not divide, the answer will be a percentage.
# for full spectrum
full_iri_R <- full_iri %>%
mutate(product = mean_grouped_refl/100 * filters.Full.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
full_noniri_male_R <- full_noniri_male %>%
mutate(product = mean_grouped_refl/100 * filters.Full.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
full_noniri_female_R <- full_noniri_female %>%
mutate(product = mean_grouped_refl/100 * filters.Full.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
# for nir
nir_iri_R <- nir_iri %>%
mutate(product = mean_grouped_refl/100 * filters.nir.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
nir_noniri_male_R <- nir_noniri_male %>%
mutate(product = mean_grouped_refl/100 * filters.nir.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
nir_noniri_female_R <- nir_noniri_female %>%
mutate(product = mean_grouped_refl/100 * filters.nir.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
# for visible (no UV)
vis_iri_R <- vis_iri %>%
mutate(product = mean_grouped_refl/100 * filters.vis.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
vis_noniri_male_R <- vis_noniri_male %>%
mutate(product = mean_grouped_refl/100 * filters.vis.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()
vis_noniri_female_R <- vis_noniri_female %>%
mutate(product = mean_grouped_refl/100 * filters.vis.f) %>% #multiply
select(wl,product) %>% # Keep only the product and wl.
as.rspec()3. Find the integral and standardize
Find the area under the curve for the product that we just calculated and the area under the curve for the raw irradiance of the illumination source. Divide the area under curve for sample by the area under curve for irradiance to get the proportion of light that each patch is able to reflect.
# calculate the reflectivity for each colour in the full spectrum
reflectivity_full_iri <- summary(full_iri_R)$B1/summary(filters[,c(1,2)])$B1
reflectivity_full_noniri_male <- summary(full_noniri_male_R)$B1/summary(filters[,c(1,2)])$B1
reflectivity_full_noniri_female <- summary(full_noniri_female_R)$B1/summary(filters[,c(1,2)])$B1
reflectivity_full <- c("iridescent" = reflectivity_full_iri, "non_iridescent_male" = reflectivity_full_noniri_male, "non_iridescent_female" = reflectivity_full_noniri_female)
# calculate the reflectivity for each colour in the nir spectrum
reflectivity_nir_iri <- summary(nir_iri_R)$B1/summary(filters[,c(1,4)])$B1
reflectivity_nir_noniri_male <- summary(nir_noniri_male_R)$B1/summary(filters[,c(1,4)])$B1
reflectivity_nir_noniri_female <- summary(nir_noniri_female_R)$B1/summary(filters[,c(1,4)])$B1
reflectivity_nir <- c("iridescent" = reflectivity_nir_iri, "non_iridescent_male" = reflectivity_nir_noniri_male, "non_iridescent_female" = reflectivity_nir_noniri_female)
# calculate the reflectivity for each colour in the visible spectrum
reflectivity_vis_iri <- summary(vis_iri_R)$B1/summary(filters[,c(1,3)])$B1
reflectivity_vis_noniri_male <- summary(vis_noniri_male_R)$B1/summary(filters[,c(1,3)])$B1
reflectivity_vis_noniri_female <- summary(vis_noniri_female_R)$B1/summary(filters[,c(1,3)])$B1
reflectivity_vis <- c("iridescent" = reflectivity_vis_iri, "non_iridescent_male" = reflectivity_vis_noniri_male, "non_iridescent_female" = reflectivity_vis_noniri_female)
# compile all reflectivirty measures, convert to data frame and save
reflectivity <- data.frame(colour = names(reflectivity_full), full = reflectivity_full, nir = reflectivity_nir, vis = reflectivity_vis)4. Visualise
Since reflectivity is a constant, we can plot it in a bar graph.
reflectivity_longer <- reflectivity %>% pivot_longer(cols = -colour, names_to = "treatment", values_to = "reflectivity") %>% arrange(treatment)
reflectivity_longer[reflectivity_longer == "non_iridescent_female"] <- "female non-iridescent"
reflectivity_longer[reflectivity_longer == "non_iridescent_male"] <- "male non-iridescent"
reflectivity_longer[reflectivity_longer == "full"] <- "Full"
reflectivity_longer[reflectivity_longer == "nir"] <- "Near-infrared"
reflectivity_longer[reflectivity_longer == "vis"] <- "Visible"
reflectivity_longer$reflectivity <- reflectivity_longer$reflectivity*100
# Reorder the factor levels of region_colour
reflectivity_longer$colour <- factor(reflectivity_longer$colour,
levels = c("female non-iridescent","male non-iridescent", "iridescent"))
ggplot(data=reflectivity_longer, aes(x=treatment, y=reflectivity, fill = colour)) +
geom_bar(position = "dodge", stat = "identity", alpha = .75) +
geom_text(aes(label = round(reflectivity, 2)), # Add labels, including rounding decimal points (round)
position = position_dodge(width = 0.9), # Adjust position to match the bars
vjust = -0.5, # Adjust vertical position (move labels above bars)
size = 3) + # Adjust the text size
coord_cartesian(ylim = c(0, 100)) +
labs(x = "Treatment", y = "Reflectivity (%)", fill = "Colour") +
scale_fill_manual(values = c("iridescent"="blue", "female non-iridescent"="orange", "male non-iridescent"="red")) +
# ylim(0,1)+
theme_bw()
Weigh reflectivity
Calculate the total reflectivity values for each bug.
# create empty columns for reflectivity values in bug datasheet
datasheet$reflectivity_full_iri <- NA
datasheet$reflectivity_full_noniri <- NA
datasheet$reflectivity_nir_iri <- NA
datasheet$reflectivity_nir_noniri <- NA
datasheet$reflectivity_vis_iri <- NA
datasheet$reflectivity_vis_noniri <- NA
# add colour category reflectivity values to bug datasheet
for(i in 1:nrow(datasheet)) {
if (datasheet[i, "sex"] == "male") {
# full
datasheet[i, "reflectivity_full_iri"] <- reflectivity[1, "full"]
datasheet[i, "reflectivity_full_noniri"] <- reflectivity[2, "full"]
# nir
datasheet[i, "reflectivity_nir_iri"] <- reflectivity[1, "nir"]
datasheet[i, "reflectivity_nir_noniri"] <- reflectivity[2, "nir"]
# visible
datasheet[i, "reflectivity_vis_iri"] <- reflectivity[1, "vis"]
datasheet[i, "reflectivity_vis_noniri"] <- reflectivity[2, "vis"]
} else {
# full
datasheet[i, "reflectivity_full_iri"] <- reflectivity[1, "full"]
datasheet[i, "reflectivity_full_noniri"] <- reflectivity[3, "full"]
# nir
datasheet[i, "reflectivity_nir_iri"] <- reflectivity[1, "nir"]
datasheet[i, "reflectivity_nir_noniri"] <- reflectivity[3, "nir"]
# visible
datasheet[i, "reflectivity_vis_iri"] <- reflectivity[1, "vis"]
datasheet[i, "reflectivity_vis_noniri"] <- reflectivity[3, "vis"]
}
}
# calculate bug reflectivity values (weighted by proportion of each colour patch - values independent of body size)
datasheet <- datasheet %>%
mutate(
full_wrfl = (proportion_iridescent*reflectivity_full_iri)+(proportion_noniridescent*reflectivity_full_noniri),
nir_wrfl = (proportion_iridescent*reflectivity_nir_iri)+(proportion_noniridescent*reflectivity_nir_noniri),
vis_wrfl = (proportion_iridescent*reflectivity_vis_iri)+(proportion_noniridescent*reflectivity_vis_noniri)
)
# calculate bug reflectivity values (weighted by area of each colour patch - values dependent on body size)
datasheet <- datasheet %>%
mutate(
full_wrfl_area = (area_iridescent_mm2*reflectivity_full_iri)+(area_noniridescent_mm2*reflectivity_full_noniri),
nir_wrfl_area = (area_iridescent_mm2*reflectivity_nir_iri)+(area_noniridescent_mm2*reflectivity_nir_noniri),
vis_wrfl_area = (area_iridescent_mm2*reflectivity_vis_iri)+(area_noniridescent_mm2*reflectivity_vis_noniri)
)Heating analysis
We conducted heating experiments using a solar simulator which approximates the spectral power distribution of natural sunlight. To understand the effect on heating, we measure two different aspects of heating:
- ΔT: the change/difference in temperature from the last point from the starting point, in which the final time point is after 5 minutes
- βT: the maximum slope of the heating curve (in degrees C per second), which is representative of the heating rate. Note that we took temperature measurements every 10 seconds, and averaged the slope over 20 seconds (2 cycles). Averaging over a longer time interval smooths out the curve to account for anomalies in heating.
Calculate heating statistics
Create a function to extract our heating summary statistics.
# function to get heating rates summary stats
ht<-function(heating,time_interval=10){ # inputs are temp of the bug at each point and the sampling time
delta<-heating[length(heating)]-heating[1] # this gives the total delta/change in bug temp
tmp_slope<-c() # create empty vector to store temp slope (change per time)
count=1
while(count<length(heating)){ # for all intervals
tmp<-(heating[count+2]-heating[count])/(time_interval*2) # get the change in temp, in degrees C per second, change count+time interval accordingly
tmp_slope<-c(tmp_slope,tmp) # save the output
count=count+1
}
slope=max(na.omit(tmp_slope)) # get the max slope (temp increase) -> note that this won't look at negative slopes
return(c(delta,slope)) # return max slope and total delta
}Run the function every bug we have data for.
# create new dataframe
reference <- datasheet
heating_data_path <- "data/heating_data/"
# add image paths for plotting
reference$imgpath <- paste0("data/bug_pics/", reference$code, ".png")
# get individual bugs
bugs <- dir(path = heating_data_path, pattern = "CSV") %>% gsub(pattern="\\w+\\_(\\w+).CSV", replacement="\\1") %>% unique()
# get treatments
treatments <- c("full", "nir", "uvvis")
# create empty columns to add delta and slope to each treatment
columnsAdd<-as.vector(outer(treatments, c("_delta","_slope"), paste0))
# add the columns to reference (the datasheet)
reference[,columnsAdd]<-0
# for loop to extract the heating summary stats
for(g in bugs){
for(t in treatments){
tmp<-read_csv(paste0(heating_data_path, t,"_",g,".CSV"),col_names=FALSE) # read in the thermodata
tmp<-tmp[1:30,] # make sure that all data is from cycle 1-30 only (or whatever cycle number)
tmp<-tmp[,3:4] # this is the bug temp and air temp
colnames(tmp)<-c("temp_bug", "temp_air")
tmp$temp_standardized <- tmp$temp_bug-tmp$temp_air # this just creates an extra column if you want to standardize the bug temp by air temp
measures<-ht(heating=tmp$temp_bug,time_interval=10) # this will give delta and slope for bug
# can control for air temp as needed, temp_bug is just bug, delta_temp_air is bug-air
pos<-grep(paste0("\\b",g,"\\b"),reference$code)
reference[pos,paste0(t,c("_delta","_slope"))]<-data.frame(measures[1],measures[2])
}
}Plot solar simulation trials
The following code was used to troubleshoot outliers and check the effect of air temperature on the solar simulator results. Here we calculate delta_temp, which standardizes the bug temperature using the starting air temperature, and delta_temp_air, which standardizes the bug temperature using the air temperature recorded simultaneously.
# empty tibble
solarsim_runs<-tibble(bug=character(),temp_bug=numeric(),temp_air=numeric(),delta_temp=numeric(),time=numeric(),trial=character(),treatment=character())
# bugs is list of individuals
for(f in bugs){
files=dir(path = heating_data_path, pattern=paste0("_",f,".CSV")) #list of file per individual ID
for(all in files){
tmp<-read_csv(paste0(heating_data_path, all),col_names=FALSE)
tmp<-tmp[1:60,]
tmp<-tmp[,c(3,4)]
colnames(tmp)<-c("temp_bug","temp_air")
tmp$delta_temp<-tmp$temp_bug-tmp$temp_bug[1] # standardizing using bug starting temp
tmp$delta_temp_air<-tmp$temp_bug-tmp$temp_air # standardizing using simultaneous air temp
tmp$time<-seq(0,590,10)
# tmp$trial<-rep(x=gsub(pattern="(\\w+)\\_\\w+\\_\\w+.CSV",replacement="\\1",x=all),times=nrow(tmp)) # get the repetition ID
tmp$treatment<-rep(x=gsub(pattern="^(\\w+)\\_(\\w+)\\.CSV",replacement="\\1",x=all),times=nrow(tmp))
tmp$bug<-f# bug ID
tmp$sex <- rep(x=gsub(pattern="^[0-9]+",replacement="",x=f),times=nrow(tmp))
solarsim_runs<-rbind(solarsim_runs,tmp)
}
}We can plot the full, NIR, and UV-visible runs from the solar simulator for each bug.
# plot the individual solar sim runs (controlled using starting air temperature)
ggplot(solarsim_runs,
aes(y=delta_temp, x=time, colour=treatment)) +
geom_line(aes()) +
# geom_line(aes(colour="black",y=temp_air-temp_air[1])) +
facet_wrap(~ bug, nrow = 5) +
# coord_cartesian(ylim = c(18, 23), xlim = c(0, 300)) +
coord_cartesian(ylim = c(-1, 20), xlim = c(0, 300)) +
# coord_cartesian(ylim = c(18, 40), xlim = c(0, 300)) +
theme_bw()
# plot the temperature of the bugs run with full spectrum light
solarsim_runs %>%
filter(
# sex == "male" | sex == "female",
treatment == "full") %>%
ggplot(aes(y = temp_bug, x = time, group = factor(bug), color = bug)) +
geom_line(alpha = 0.5, linewidth = 1) +
theme_bw()
Statistical analysis
Size-dependent reflectivity and heating
To check that our setup could detect variation in heating rates, we used a linear regression to examine the relationship between heating and area-dependent reflectivity, which accounts for the reflectivity of absolute surface area of the bug. We calculated area-dependent reflectivity as the reflectivity of the colour patches multiplied by their area – as opposed to the proportion – which means this is a measure of reflectivity that incorporates body size.
Model summaries
Results of the linear models.
models_sizerefl <- list(
"full" = lm(data = reference, full_delta ~ full_wrfl_area),
"NIR" = lm(data = reference, nir_delta ~ nir_wrfl_area),
"UV-VIS" = lm(data = reference, uvvis_delta ~ vis_wrfl_area),
"full" = lm(data = reference, full_slope ~ full_wrfl_area),
"NIR" = lm(data = reference, nir_slope ~ nir_wrfl_area),
"UV-VIS" = lm(data = reference, uvvis_slope ~vis_wrfl_area)
)
modelsummary(models_sizerefl, stars = TRUE,
estimate = "{estimate}",
statistic = c("{std.error}",
"{statistic}",
"{p.value}"),
output = "tinytable") %>%
group_tt(j = list("Delta (ΔT)" = 2:4, "Slope (βT)" = 5:7)) #%>% save_tt("models_sizerefl.docx")| Delta (ΔT) | Slope (βT) | |||||
|---|---|---|---|---|---|---|
| full | NIR | UV-VIS | full | NIR | UV-VIS | |
| + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001 | ||||||
| (Intercept) | 15.093 | 8.958 | 5.025 | 0.169 | 0.095 | 0.056 |
| 0.892 | 0.487 | 0.336 | 0.015 | 0.008 | 0.005 | |
| 16.929 | 18.403 | 14.966 | 11.044 | 11.482 | 10.539 | |
| <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | |
| full_wrfl_area | -0.047 | -0.001 | ||||
| 0.028 | 0.000 | |||||
| -1.697 | -2.636 | |||||
| 0.097 | 0.012 | |||||
| nir_wrfl_area | -0.036 | -0.001 | ||||
| 0.012 | 0.000 | |||||
| -3.126 | -3.091 | |||||
| 0.003 | 0.004 | |||||
| vis_wrfl_area | -0.029 | -0.001 | ||||
| 0.018 | 0.000 | |||||
| -1.610 | -1.999 | |||||
| 0.115 | 0.052 | |||||
| Num.Obs. | 43 | 43 | 43 | 43 | 43 | 43 |
| R2 | 0.066 | 0.193 | 0.059 | 0.145 | 0.189 | 0.089 |
| R2 Adj. | 0.043 | 0.173 | 0.036 | 0.124 | 0.169 | 0.067 |
| AIC | 201.4 | 149.8 | 115.8 | -148.4 | -200.3 | -241.0 |
| BIC | 206.6 | 155.1 | 121.0 | -143.1 | -195.0 | -235.7 |
| Log.Lik. | -97.675 | -71.898 | -54.881 | 77.200 | 103.151 | 123.509 |
| RMSE | 2.35 | 1.29 | 0.87 | 0.04 | 0.02 | 0.01 |
Visualise
Models with ΔT (final temperature after 5 minutes) as response variable.
# treatment = full
ggplot(data = reference, aes(y=full_delta, x=full_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = nir
ggplot(data = reference, aes(y=nir_delta, x=nir_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = uvvis
ggplot(data = reference, aes(y=uvvis_delta, x=vis_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
Models with βT (heating rate) as response variable.
# models with slope as response
# treatment = full
ggplot(data = reference, aes(y=full_slope, x=full_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = nir
ggplot(data = reference, aes(y=nir_slope, x=nir_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = uvvis
ggplot(data = reference, aes(y=uvvis_slope, x=vis_wrfl_area)) +
geom_smooth(method = lm) +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "Reflectivity (%)",
colour = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
Body size and reflectivity
We examined the relationship between bug reflectivity (mean reflectivity of each colour patch weighted by the proportion of that colour) and body size, our two predictor variables.
Model summaries
Check the relationship between our two predictor variables: reflectivity and body size (measured as the area in mm2 of the dorsal surface) and check the correlation coefficient. Note that body size was square root transformed.
# correlation plot
# subset all the columns you need
corr_plot <- subset(reference, select = c(
"full_wrfl",
"nir_wrfl",
"vis_wrfl",
"area_mm2"
)
) %>% na.omit() # remove NA values
M = cor(corr_plot, method = "spearman")
corrplot(M, method = 'number') # colourful numbers
Visualise
# linear relationship between body size and colour (Figure 4a)
lm(data = reference, full_wrfl ~ sqrt(area_mm2)) %>% summary() # summary of model##
## Call:
## lm(formula = full_wrfl ~ sqrt(area_mm2), data = reference)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.077007 -0.030367 -0.003526 0.020152 0.114974
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.076185 0.054591 -1.396 0.17
## sqrt(area_mm2) 0.025481 0.004757 5.356 3.54e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04354 on 41 degrees of freedom
## Multiple R-squared: 0.4117, Adjusted R-squared: 0.3973
## F-statistic: 28.69 on 1 and 41 DF, p-value: 3.545e-06# lm(data = reference, full_wrfl ~ sqrt(area_mm2)) %>% plot() # diagnostic plot
ggplot(data = reference, aes(x=sqrt(area_mm2), y=full_wrfl)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "Reflectivity (%)", x = "Body size (mm²)") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# male vs female size
ggplot(data = reference, aes(x=sex, y=area_mm2)) +
geom_boxplot() +
# geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "Body size (mm²)", x = "Sex") +
# scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
# scale_colour_manual(values = c("male" = "#7570b3", "female" = "#1b9e77")) +
theme_bw()
# male vs female reflectivity
ggplot(data = reference, aes(x=sex, y=full_wrfl*100)) +
geom_boxplot() +
# geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "Reflectivity (%)", x = "Sex") +
# scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
# scale_colour_manual(values = c("male" = "#7570b3", "female" = "#1b9e77")) +
theme_bw()
Check the variation in reflectivity between sexes, morphs, and in different wavelengths.
wrfl_variation <- reference %>% select(code, sex, proportion_iridescent, full_wrfl, nir_wrfl, vis_wrfl, imgpath)
wrfl_variation <- wrfl_variation %>%
mutate(morph = ifelse(sex == "male" & proportion_iridescent <= .50, "male (red-orange morph)",
ifelse(sex == "male" & proportion_iridescent > .50, "male (blue morph)",
ifelse(sex == "female", "female", NA))))
wrfl_variation <- wrfl_variation %>% pivot_longer(cols = c(full_wrfl, nir_wrfl, vis_wrfl),
names_to = "spectra", names_pattern = "(.*)_wrfl",
values_to = "wrfl")
wrfl_variation$morph <- factor(wrfl_variation$morph, levels = c("female", "male (red-orange morph)", "male (blue morph)")) # Specify the desired order
wrfl_variation$spectra <- factor(wrfl_variation$spectra, levels = c("vis", "nir", "full")) # Specify the desired order
# plot
ggplot(data = wrfl_variation, aes(x=spectra, y=wrfl*100)) +
geom_boxplot() +
geom_point(aes(color = morph), alpha = .70, size = 3,
position = position_jitter(seed = 23)) +
# geom_image(aes(image=imgpath), size = .07) +
scale_x_discrete(labels = c("full" = "full", "nir" = "near-infrared", "vis" = "visible")) +
labs(y = "Reflectivity (%)", x = "Wavelength range", colour = "Sex") +
scale_colour_manual(values = c("male (blue morph)" = "#2A788EFF", "male (red-orange morph)" = "#440154FF", "female" = "#3dc18c")) +
theme_bw()
Size-independent reflectivity and heating
We can run a regression model to address the collinearity between predictors by using the residuals of one predictor regressed on the other. To to this we can: 1. Fit a regression model with one predictor (e.g., color) regressed on the other (e.g., size). 2. Extract the residuals from this regression. 3. Use these residuals as the new predictor in your final regression model with the response variable (heating rate).
Calculate colour residuals.
# add residuals to reference data frame
reference$full_color_residuals <- lm(data = reference, full_wrfl ~ sqrt(area_mm2)) %>% residuals()
reference$nir_color_residuals <- lm(data = reference, nir_wrfl ~ sqrt(area_mm2)) %>% residuals()
reference$vis_color_residuals <- lm(data = reference, vis_wrfl ~ sqrt(area_mm2)) %>% residuals()
# check to see that residuals are now independent of body size
ggplot(data = reference, aes(y=full_color_residuals, x=sqrt(area_mm2))) +
geom_smooth(method = lm) +
geom_point(aes(color = sex)) +
labs(x = "Reflectivity (controlled by size)", y = "Square root(Body size in mm²)",
color = "Sex") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
### Model summaries
Results of the linear models.
# create the models
models <- list(
"full" = lm(data = reference, full_delta ~ full_color_residuals + sqrt(area_mm2)),
"NIR" = lm(data = reference, nir_delta ~ nir_color_residuals + sqrt(area_mm2)),
"UV-VIS" = lm(data = reference, uvvis_delta ~ vis_color_residuals + sqrt(area_mm2)),
"full" = lm(data = reference, full_slope ~ full_color_residuals + sqrt(area_mm2)),
"NIR" = lm(data = reference, nir_slope ~ nir_color_residuals + sqrt(area_mm2)),
"UV-VIS" = lm(data = reference, uvvis_slope ~ vis_color_residuals + sqrt(area_mm2))
)
modelsummary(models, stars = TRUE,
estimate = "{estimate}",
statistic = c("{std.error}",
"{statistic}",
"{p.value}"),
output = "tinytable") %>%
group_tt(j = list("Delta (ΔT)" = 2:4, "Slope (βT)" = 5:7)) #%>% save_tt("modeltable_UPDATED.docx")| Delta (ΔT) | Slope (βT) | |||||
|---|---|---|---|---|---|---|
| full | NIR | UV-VIS | full | NIR | UV-VIS | |
| + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001 | ||||||
| (Intercept) | 19.098 | 11.828 | 6.115 | 0.270 | 0.149 | 0.076 |
| 3.023 | 1.666 | 1.128 | 0.052 | 0.029 | 0.018 | |
| 6.318 | 7.100 | 5.421 | 5.226 | 5.209 | 4.277 | |
| <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | |
| full_color_residuals | -5.073 | -0.102 | ||||
| 8.647 | 0.148 | |||||
| -0.587 | -0.691 | |||||
| 0.561 | 0.494 | |||||
| sqrt(area_mm2) | -0.473 | -0.374 | -0.139 | -0.012 | -0.007 | -0.003 |
| 0.263 | 0.145 | 0.098 | 0.005 | 0.002 | 0.002 | |
| -1.795 | -2.574 | -1.415 | -2.694 | -2.716 | -1.686 | |
| 0.080 | 0.014 | 0.165 | 0.010 | 0.010 | 0.100 | |
| nir_color_residuals | -6.384 | -0.077 | ||||
| 3.458 | 0.059 | |||||
| -1.846 | -1.299 | |||||
| 0.072 | 0.201 | |||||
| vis_color_residuals | -3.539 | -0.095 | ||||
| 5.160 | 0.081 | |||||
| -0.686 | -1.175 | |||||
| 0.497 | 0.247 | |||||
| Num.Obs. | 43 | 43 | 43 | 43 | 43 | 43 |
| R2 | 0.082 | 0.201 | 0.058 | 0.162 | 0.185 | 0.095 |
| R2 Adj. | 0.036 | 0.161 | 0.011 | 0.120 | 0.144 | 0.050 |
| AIC | 202.6 | 151.4 | 117.8 | -147.3 | -198.1 | -239.3 |
| BIC | 209.6 | 158.4 | 124.9 | -140.2 | -191.0 | -232.3 |
| Log.Lik. | -97.299 | -71.684 | -54.908 | 77.633 | 103.039 | 123.666 |
| RMSE | 2.33 | 1.28 | 0.87 | 0.04 | 0.02 | 0.01 |
Visualise
Models with ΔT (final temperature after 5 minutes) as response variable.
# treatment = full
ggplot(data = reference, aes(y=full_delta, x=full_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = nir
ggplot(data = reference, aes(y=nir_delta, x=nir_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = uvvis
ggplot(data = reference, aes(y=uvvis_delta, x=vis_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "ΔT (°C)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
Models with βT (heating rate) as response variable.
# treatment = full
ggplot(data = reference, aes(y=full_slope, x=full_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = nir
ggplot(data = reference, aes(y=nir_slope, x=nir_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
# treatment = uvvis
ggplot(data = reference, aes(y=uvvis_slope, x=vis_color_residuals)) +
geom_smooth(method = lm, color = "black") +
geom_point(aes(color = sex), alpha = .25, size = 7) +
geom_image(aes(image=imgpath), size = .07) +
labs(y = "βT (Δ°C/second)", x = "more iridescent Reflectivity (controlled by size) less iridescent") +
scale_colour_manual(values = c("male" = "#571267", "female" = "#3dc18c")) +
theme_bw()## `geom_smooth()` using formula = 'y ~ x'
Separate models
Body size only
models_size <- list(
"full" = lm(data = reference, full_delta ~ sqrt(area_mm2)),
"NIR" = lm(data = reference, nir_delta ~ sqrt(area_mm2)),
"UV-VIS" = lm(data = reference, uvvis_delta ~ sqrt(area_mm2)),
"full" = lm(data = reference, full_slope ~ sqrt(area_mm2)),
"NIR" = lm(data = reference, nir_slope ~ sqrt(area_mm2)),
"UV-VIS" = lm(data = reference, uvvis_slope ~ sqrt(area_mm2))
)
modelsummary(models_size, stars = TRUE,
estimate = "{estimate}",
statistic = c("{std.error}",
"{statistic}",
"{p.value}"),
output = "tinytable") %>%
group_tt(j = list("Delta (ΔT)" = 2:4, "Slope (βT)" = 5:7)) | Delta (ΔT) | Slope (βT) | |||||
|---|---|---|---|---|---|---|
| full | NIR | UV-VIS | full | NIR | UV-VIS | |
| + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001 | ||||||
| (Intercept) | 19.098 | 11.828 | 6.115 | 0.270 | 0.149 | 0.076 |
| 2.998 | 1.714 | 1.121 | 0.051 | 0.029 | 0.018 | |
| 6.370 | 6.900 | 5.457 | 5.259 | 5.166 | 4.257 | |
| <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | |
| sqrt(area_mm2) | -0.473 | -0.374 | -0.139 | -0.012 | -0.007 | -0.003 |
| 0.261 | 0.149 | 0.098 | 0.004 | 0.003 | 0.002 | |
| -1.809 | -2.501 | -1.425 | -2.711 | -2.694 | -1.678 | |
| 0.078 | 0.016 | 0.162 | 0.010 | 0.010 | 0.101 | |
| Num.Obs. | 43 | 43 | 43 | 43 | 43 | 43 |
| R2 | 0.074 | 0.132 | 0.047 | 0.152 | 0.150 | 0.064 |
| R2 Adj. | 0.051 | 0.111 | 0.024 | 0.131 | 0.130 | 0.041 |
| AIC | 201.0 | 152.9 | 116.3 | -148.8 | -198.3 | -239.9 |
| BIC | 206.2 | 158.2 | 121.6 | -143.5 | -193.0 | -234.6 |
| Log.Lik. | -97.483 | -73.442 | -55.160 | 77.378 | 102.151 | 122.937 |
| RMSE | 2.34 | 1.34 | 0.87 | 0.04 | 0.02 | 0.01 |
Residuals only
# ONLY RESIDUALS
models_residuals <- list(
"full" = lm(data = reference, full_delta ~ full_color_residuals),
"NIR" = lm(data = reference, nir_delta ~ nir_color_residuals),
"UV-VIS" = lm(data = reference, uvvis_delta ~ vis_color_residuals),
"full" = lm(data = reference, full_slope ~ full_color_residuals),
"NIR" = lm(data = reference, nir_slope ~ nir_color_residuals),
"UV-VIS" = lm(data = reference, uvvis_slope ~ vis_color_residuals)
)
modelsummary(models_residuals, stars = TRUE,
estimate = "{estimate}",
statistic = c("{std.error}",
"{statistic}",
"{p.value}"),
output = "tinytable")| full | NIR | UV-VIS | full | NIR | UV-VIS | |
|---|---|---|---|---|---|---|
| + p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001 | ||||||
| (Intercept) | 13.714 | 7.572 | 4.530 | 0.132 | 0.072 | 0.046 |
| 0.377 | 0.216 | 0.139 | 0.007 | 0.004 | 0.002 | |
| 36.330 | 35.041 | 32.627 | 19.544 | 19.219 | 20.943 | |
| <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | <0.001 | |
| full_color_residuals | -5.073 | -0.102 | ||||
| 8.878 | 0.159 | |||||
| -0.571 | -0.644 | |||||
| 0.571 | 0.523 | |||||
| nir_color_residuals | -6.384 | -0.077 | ||||
| 3.687 | 0.064 | |||||
| -1.731 | -1.208 | |||||
| 0.091 | 0.234 | |||||
| vis_color_residuals | -3.539 | -0.095 | ||||
| 5.223 | 0.083 | |||||
| -0.678 | -1.150 | |||||
| 0.502 | 0.257 | |||||
| Num.Obs. | 43 | 43 | 43 | 43 | 43 | 43 |
| R2 | 0.008 | 0.068 | 0.011 | 0.010 | 0.034 | 0.031 |
| R2 Adj. | -0.016 | 0.045 | -0.013 | -0.014 | 0.011 | 0.008 |
| AIC | 203.9 | 156.0 | 117.9 | -142.1 | -192.8 | -238.4 |
| BIC | 209.2 | 161.2 | 123.2 | -136.8 | -187.5 | -233.1 |
| Log.Lik. | -98.964 | -74.979 | -55.959 | 74.049 | 99.399 | 122.191 |
| RMSE | 2.42 | 1.38 | 0.89 | 0.04 | 0.02 | 0.01 |